\section{轨道摄动}

\begin{frame}{\thesection.\ \insertsection}
开普勒二体轨道是高度理想化的运动模型。
\begin{itemize}
    \item 短期预测效果良好，但实际运动受多种因素影响会偏离开普勒轨道。
\end{itemize}
\begin{center}\includegraphics{fig_5_p53.pdf}\end{center}
\end{frame}

\begin{frame}{\thesection.\ \insertsection}
主要摄动因素包括：
\begin{enumerate}
    \item 主天体 \( m_1 \) 的非球形（地球并非完美球体）
    \begin{itemize}
        \item 开普勒轨道假设 \( m_1 \) 和 \( m_2 \) 为质点
        \begin{itemize}
            \item 若主天体 \( m_1 \) 具有球对称性，则可视为质点
            \item 由于 \( m_2 \) 尺寸远小于其与 \( m_1 \) 中心的距离，此近似对 \( m_2 \) 成立
        \end{itemize}
    \end{itemize}
    \item 其他天体的引力场影响
    \begin{itemize}
        \item 地球轨道航天器还受太阳、月球引力影响（其他行星影响更弱）
    \end{itemize}
    \item 大气阻力
    \begin{itemize}
        \item 近地轨道存在残余大气，会对航天器产生阻力导致轨道衰减
    \end{itemize}
    \item 太阳光压
    \begin{itemize}
        \item 太阳光子撞击航天器受光面产生压力，源于光子与表面的动量交换
    \end{itemize}
\end{enumerate}
\end{frame}

\begin{frame}{\thesection.\ \insertsection}
这些效应构成对开普勒轨道的摄动。\\
设 $\vec f_p$ 为摄动加速度，则真实运动方程为：
\[ \ddot{\vec{r}} = -\frac{\mu}{r^3} \vec{r} + \vec{f}_p \]
处理摄动有两种方法：
\begin{enumerate}
    \item 特殊摄动法
    \begin{itemize}
        \item 通过数值积分计算摄动影响
        \item 解仅适用于特定初始条件
    \end{itemize}
    \item 普遍摄动法
    \begin{itemize}
        \item 解析求解摄动影响，推导轨道要素变化率：
        \[\frac{da}{dt}, \frac{de}{dt}, \frac{di}{dt}, \frac{d\Omega}{dt}, \frac{d\omega}{dt}\]
        \item 推导过程常需近似处理，精度低于特殊摄动法
    \end{itemize}
\end{enumerate}
\end{frame}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsection{特殊摄动法}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\begin{frame}{\thesection.\ \insertsection \\ \small\thesection.\thesubsection\ \insertsubsection}
1. 考威尔方法（20世纪初）
\begin{itemize}
    \item 考威尔方法直接数值积分运动方程
    \item 在惯性参考系 \( \mathcal{F}_I \) 中表示所有量：
    \[\vec{r} = \vec{F}^T_I r, \quad \vec{v} = \vec{F}^T_I v, \quad \vec{f}_p = \vec{F}^T_I f_p, \quad \vec{r}_0 = \vec{F}^T_I r_0, \quad \vec{v}_0 = \vec{F}^T_I v_0\]
    \item \textcolor{blue}{可将运动方程改写为适合数值积分的一阶形式：}
    \[\begin{bmatrix}
    \dot{r} \\
    \dot{v}
    \end{bmatrix} =
    \begin{bmatrix}
    -\frac{\mu}{(r^\text Tr)^\frac{3}{2}}r + f_p \\
    \end{bmatrix},
    \begin{bmatrix}
    r(0) \\
    v(0)
    \end{bmatrix} =
    \begin{bmatrix}
    r_0 \\
    v_0
    \end{bmatrix}\]
\end{itemize}
注记：
\begin{itemize}
    \item 为保证精度通常需要很小的时间步长，计算成本高
    \item 数值积分的舍入误差会快速累积，导致长期预测不准确
    \item 随着计算能力和数值技术的发展，该问题已较过去有所改善
\end{itemize}
\end{frame}

\begin{frame}{\thesection.\ \insertsection \\ \small\thesection.\thesubsection\ \insertsubsection}
2. 恩克方法（1857年）
该方法通过数值积分真实（受摄）轨道与参考开普勒轨道的偏差来求解。真实轨道满足：
\[\ddot{\vec{r}} = -\frac{\mu}{r^3} \vec{r} + \vec{f}_p, \quad \vec{r}(0) = \vec{r}_0, \quad \dot{\vec{r}}(0) = \vec{v}_0\]
参考开普勒轨道为：
\[\ddot{\vec{p}} = -\frac{\mu}{\rho^3} \vec{\rho}, \quad \vec{\rho}(0) = \vec{r}_0, \quad \dot{\vec{\rho}}(0) = \vec{v}_0\]
定义轨道偏差为：
\[\delta \vec{r} \triangleq \vec{r} - \vec{\rho}\]
\textcolor{blue}{两式相减得到偏差方程：}
\[\delta \ddot{\vec{r}} = -\frac{\mu}{\rho^3} \left[ \delta \vec{r} - \left( 1 - \frac{\rho^3}{r^3} \right) \vec{r} \right] + \vec{f}_p, \quad \delta \vec{r}(0) = \vec{0}, \quad \delta \dot{\vec{r}}(0) = \vec{0}\]
\[\vec{r} = \vec{\rho} + \delta \vec{r}\]
\end{frame}

\begin{frame}{\thesection.\ \insertsection \\ \small\thesection.\thesubsection\ \insertsubsection}
注记：
\begin{itemize}
    \item 该方法计算效率可达考威尔方法的10倍
\begin{itemize}
    \item 由于 \(\delta \vec{r}\) 的变化比 \(\vec{r}\) 缓慢得多，在相同精度要求下可使用更大的时间步长
\end{itemize}
\end{itemize}
\end{frame}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsection{普遍摄动法}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\begin{frame}{\thesection.\ \insertsection \\ \small\thesection.\thesubsection\ \insertsubsection}
该方法推导轨道要素变化的解析表达式。\\
采用柱坐标系 \( \mathcal{F}_o \)，其基向量为 \( \vec{x}_o, \vec{y}_o \) 和 \( \vec{z}_o \)。
\begin{center}\includegraphics{fig_5_1.pdf}\end{center}
\begin{center}\textcolor{blue}{图 \arabic{section}.1:} 柱坐标系\end{center}
将摄动力在 \( \mathcal{F}_o \) 坐标系中表示为：
\[\vec{f}_p = f_r \vec{x}_o + f_\theta \vec{y}_o + f_z \vec{z}_o\]
\end{frame}

\begin{frame}{\thesection.\ \insertsection \\ \small\thesection.\thesubsection\ \insertsubsection}
可得到轨道要素的变分方程：
\begin{align*}
\frac{da}{dt} &= \frac{2a^2}{\sqrt{\mu a(1 - e^2)}} [e \sin \theta f_r + (1 + e \cos \theta) f_\theta] \\
\frac{de}{dt} &= \sqrt{\frac{a(1 - e^2)}{\mu}} \left[\sin \theta f_r + \frac{2 \cos \theta + e(1 + \cos^2 \theta)}{1 + e \cos \theta} f_\theta\right] \\
\frac{di}{dt} &= \sqrt{\frac{a(1 - e^2)}{\mu}} \frac{\cos(\omega + \theta)}{1 + e \cos \theta} f_z \\
\frac{d\Omega}{dt} &= \sqrt{\frac{a(1 - e^2)}{\mu}} \frac{\sin(\omega + \theta)}{\sin i(1 + e \cos \theta)} f_z \\
\frac{d\omega}{dt} &= \sqrt{\frac{a(1 - e^2)}{\mu}} \left[ -\frac{\cos \theta}{e} + \frac{(2 + e \cos \theta) \sin \theta}{e(1 + e \cos \theta)} f_\theta - \frac{\sin(\omega + \theta)}{\tan i(1 + e \cos \theta)} f_z \right]
\end{align*}
\end{frame}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsection{非球形主天体的引力摄动}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\begin{frame}{\thesection.\ \insertsection \\ \small\thesection.\thesubsection\ \insertsubsection}
建立摄动方法后，需要给出摄动加速度 \( \vec{f}_p \) 的具体表达式。
\begin{itemize}
    \item 我们将重点讨论非球形主天体的影响，这对地球轨道卫星尤为重要。
    \item \textcolor{blue}{地球呈扁球体形状，两极稍扁。}
\end{itemize}
\begin{center}\includegraphics{fig_5_2.pdf}\end{center}
\begin{center}\textcolor{blue}{图 \arabic{section}.2:} 地球两极扁平示意图\end{center}
\end{frame}

\begin{frame}{\thesection.\ \insertsection \\ \small\thesection.\thesubsection\ \insertsubsection}
我们可以用柱坐标系表示单位质量的摄动力：
\[ \vec{f}_p = f_r \vec{x}_o + f_\theta \vec{y}_o + f_z \vec{z}_o \]
其中
\begin{align*}
f_r &= \frac{3\mu J_2 R_e^2}{2r^4} (3\sin^2 i \sin^2 (\omega + \theta) - 1) \\
f_\theta &= \frac{3\mu J_2 R_e^2}{2r^4} \sin^2 i \sin (2(\omega + \theta)) \\
f_z &= \frac{3\mu J_2 R_e^2}{2r^4} \sin 2i \sin (\omega + \theta) \\
J_2 &= 1.083 \times 10^{-3}
\end{align*}
\end{frame}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsection{$J_2$项对轨道要素的影响}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\begin{frame}{\thesection.\ \insertsection \\ \small\thesection.\thesubsection\ \insertsubsection}
得到$J_2$摄动的单位质量力表达式后，可用普遍摄动法分析其影响。
\begin{itemize}
    \item 可求得轨道要素的平均变化率：
\end{itemize}
\begin{align*}
&< \dot{\Omega} > = -\frac{3J_2R_e^2}{2(1 - e^2)^2} \sqrt{\frac{\mu}{a^7}} \cos i \\
&< \dot{\omega} > = \frac{3J_2R_e^2}{4(1 - e^2)^2} \sqrt{\frac{\mu}{a^7}} (5 \cos^2 i - 1) \\
&< \dot{a} > = < \dot{e} > = < \dot{i} > = 0
\end{align*}
\end{frame}

\begin{frame}{\thesection.\ \insertsection \\ \small\thesection.\thesubsection\ \insertsubsection}
地球扁率效应（通过$J_2$项）长期仅影响$\Omega$和$\omega$。
\begin{itemize}
    \item 轨道平面以平均速率$\langle \dot{\Omega} \rangle$绕地球自转轴旋转
    \item 近地点幅角以平均速率$\langle \dot{\omega} \rangle$绕轨道法向旋转
\end{itemize}
\begin{center}\includegraphics{fig_2_8.pdf}\end{center}
\begin{center}\textcolor{blue}{图 \arabic{section}.3:} 轨道要素示意图\end{center}
\end{frame}
